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Abstract 

The chemotaxis pathway of the bacterium Rhodobacter sphaeroides has many sim- 
ilarities to the well-studied pathway in Escherichia coli. It exhibits robust adaptation 
and has several homologues of the latter's chemotaxis proteins. Recent theoretical re- 
sults have been able to correctly predict that the chemotactic response of Escherichia 
coli exhibits the same output behavior in response to scaled ligand inputs, a dynamic 
property known as fold-change detection (FCD), or input-scale invariance. In this 
paper, we present theoretical assumptions on the R. sphaeroides chemotaxis sensing 
dynamics that can be analytically shown to yield FCD behavior in a specific ligand 
concentration range. Based on these assumptions, we construct two models of the full 
chemotaxis pathway that are able to reproduce experimental time-series data from ear- 
lier studies. To test the validity of our assumptions, we propose a series of experiments 
in which our models predict robust FCD behavior where earlier models do not. In 
this way, we illustrate how a dynamic phenotype such as FCD can be used for the 
purposes of discriminating between two models that reproduce the same experimental 
time-series data. 

1 Introduction 

Dynamic models of biological mechanisms are meaningful if they can explain experimental 
data, make a priori predictions of biological behavior and be liable to invalidation through 
testing. 

Although several competing models of a given mechanism can often be made to repro- 
duce experimental data through sufficient parameterization and tuning, in many cases it 
is possible to discriminate between such models by comparing the experimentally observed 
output response and the simulated response to a judiciously designed perturbation. This 
paper is a study of the use of a particular dynamic phenotype for the purposes of model 
discrimination. Dynamic phenotypes are distinctive, qualitative, dynamic output responses 
that are robustly maintained under a range of experimental conditions. 

An example of a dynamic phenotype is adaptation, where a system initially at steady- 
state reacts to input stimuli and then restores its pre-stimulus equilibrium. It has been shown 
that integral control is the structural feature responsible for this behavior [15]. Weber's 
law, whereby a system exhibits the same maximal amplitude in its response to two different 
inputs that are positive linear scalings of each other [5] is another example of a dynamic 
phenotype. 

This paper deals with a third dynamic phenotype, termed fold change detection (FCD) 
[5], or scale invariance. A system is said to exhibit FCD if its output responses to two 
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different input stimuli that are positive linear scalings of each other are identical (which 
makes this a stronger property than Weber's law). 

In a study [S1[T], it was predicted that the chcmotaxis system of Escherichia coli, modeled 
in |10j . would exhibit the FCD property, and these predictions were later confirmed as 
accurate [2]. The key assumption of this model, which leads to FCD, is the allosteric 
signaling structure of the methyl-accepting chcmotaxis protein receptors. 

Although significantly more complex, the chcmotaxis system of the bacterium Rhodobac- 
ter sphaeroides has many similarities to that of E. coli. It features two, rather than one, 
sensory clusters; one at the cell membrane and the other in the cytoplasm. Whilst the 
membrane cluster, as in E. coli, detects external ligand, it is as yet unknown exactly what 
the cytoplasmic cluster senses [3]. Besides detecting internalized ligand concentrations, it 
may also sense internal signals, such as signals reporting the cell's metabolic state. This 
bacterium also has multiple homologues of the E. coli chemotaxis proteins, which play 
roles similar to those found in the latter, although the exact structure of their connec- 
tivity with the two sensory clusters and the flagellum is not known with certainty. The 
ChcA homologues transduce the receptor activity to the other chemotaxis proteins through 
phosphotransfer, the CheR and CheB homologues respectively mcthylate and dcmcthylatc 
receptors, whilst the CheY proteins are believed to have a role in varying the stopping 
frequency of the bacterium's single flagellum [3]. 

Recent studies have used a model invalidation technique to suggest possible connectivities 
for the CheY proteins [5J and the CheB proteins |2J. However, upon simulation it becomes 
evident that these models do not exhibit the FCD behavior observed in E. coli. This 
suggests the question: given the similarities between the two chemotaxis pathways, does the 
R. sphaeroides chemotaxis response show FCD as does that of E. colli 

In this paper, we model the dynamics of the two R. sphaeroides receptor clusters using 
the MWC allosteric model [E] that has been used to model the receptor activity in E. coli 
in [T3J HH [TO] . We present a theorem that shows that if this is an accurate model of the 
receptor dynamics, then the receptor activities will exhibit FCD. What is more, this observed 
behavior is robust to the connectivity between the chcmotaxis proteins, the receptors and the 
flagellum. To illustrate this point, we construct two models of the integrated R. sphaeroides 
chemotaxis pathway based on our receptor dynamics assumptions, with each model featuring 
a different connectivity. We show that, in addition to reproducing previously published 
experimental data, these models also display FCD in their flagellar responses in certain 
ligand concentration ranges. Since flagellar outputs can be easily measured using tethered 
cell assays, we then suggest a series of experiments that can be used to test whether the 
models we present here are accurate compared to previously published models based on 
whether or not the flagellar response exhibits FCD. 

This work therefore makes the case that qualitative dynamic behavior could be a powerful 
property to test when discriminating between competing models. A systematic way of model 
discrimination using this approach would start with the construction of a dynamic model 
that explains experimental data. The next step would be to use the model to mathematically 
identify experimentally implementable conditions under which the system can be expected to 
exhibit a certain dynamic phenotype. The final step would be to experimentally implement 
those conditions and to compare the measured results against what is predicted in silico. 
In this way, two models which explain experimental data equally well can be discriminated 
using their dynamic phenotypes. 
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1.1 Background 



We can decompose the R. sphaeroides chemotaxis pathway into three modules, as illustrated 
in Figure [T] The sensing module includes two receptor clusters. One of these resides at the 
cell membrane and senses the concentration of external ligands L, as illustrated in Figure 
[5J The other cluster resides within the cytoplasm and measures an internalized ligand 
concentration L. Henceforth, the ~ notation will be used to denote signals associated with 
the cytoplasmic cluster. 

Receptor activity regulation 
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Figure 1: Schematic of the R. sphaeroides chemotaxis pathway. 



The dynamics of the two receptor clusters are modeled as two first-order systems. The 
membrane receptor cluster is assumed to have state to (its receptor methylation level) and 
output a (the receptor activity level). Similarly the cytoplasmic cluster has methylation 
level to as its state and its activity level a as its output. The state-space representation of 
this system is then 

m = F(a, w) 
a =G(m,L) 
m = F(a, w) 
a = G(m, L) 

where w, w are functions of the concentrations of the phorphorylatcd chemotaxis proteins 
within the cell. These functions represent the interactions between the internal state of 
the cell and the receptors. For example, w and w can represent the demcthylation of the 
receptors by the proteins CheBi,CheB2 or their methylation by the proteins ChcR,2,ChcR3. 

The cytoplasmic cluster is believed to integrate the extra-cellular ligand concentration 
L with internal cell signals. We represent these internal cell signals by u, a function of 
the concentrations of the phorphorylated chemotaxis proteins. The signal L in Figure [2] is 
assumed to have the following relation with the externally sensed ligand. 

Assumption 1. The internalized ligand concentration L is related to the external ligand 
concentration L through a linear, time invariant filter 

i = A£ + B{u)L 1 ', f el" 
L = C£, + D{u)L v 

where A g M. nxn , B : K ->• R n , C g R lxn , D : M -> E and v g R. 

With Assumption [TJ the internalized ligand concentration L can represent a variety of 
signals, including, for example, a static map that combines the externally sensed ligands L 
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Figure 2: The R. sphaeroides signalling network. 



with the internal chemotaxis protein signals u, or it can be a phase-delayed version of L or 
even, to allow for a degree of possible cooperativity, a power of L. 

In the transduction sub-system of Figure [TJ auto-phosphorylation of the chemotaxis pro- 
tein CheA2 is accelerated by the membrane cluster activity, whilst that of the CheAaA4 
complex is catalyzed by cytoplasmic cluster activity (as shown in Figure [2]) . The proteins 
CheY3, CheY4, CheYg, CheBi and ChcB2 all compete for phosphoryl groups from CheA2, 
whilst CI1CB2 and ChcYg do so from CheAaA^ The reaction rates for all of these phos- 
phorylations are given in [TH [5] . We represent this phosphotransfer network as a general 
nonlinear system, with state vector 

x=[A 2p Y 3p Y 4p (A 3 A 4 ) P Y 6p B lp B 2p f 

the individual states being the concentrations of the phosphorylated chemotaxis proteins. 
The transduction system takes as its inputs the receptor activities a, a: 

x = H(x,a,a) (2) 

where H(x,a, a) is given by the ODEs (9)-(15) in jS]. 

The outputs of this system are signals u>(x), iD(x), u(x, a), which feed back into the 
sensing subsystem, as described above. The interconnection of the phosphotransfer network 
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@ with the receptor dynamics is illustrated in Figure [H and the interconnection between 
the two subsystems can thus be written as 

m = F(a, u>(x)), a = G(m, L) 

fh = F(a, w(x)), a = G(rh,L) . . 

i = At + B(u(x, a))L v , L = C£ + D{u)L v 1 j 

x = if(x, a, a) 



External ligand L 



lu (x) 
m(x, a) 



Receptor methylation and activity 



rh = F(a, w), a = 


G(m, L) 


m = F(a, w), a = 
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L = C( + D(u)V 
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Figure 3: The interconnection of the receptors' sensing dynamics with the phosphotransfer 
network. 



As shown in Figure [5J the protein CheYg-P, possibly acting together with one or both 
of CheY3-P and CheY4-P, is believed to bind with the flagellar motor proteins to inhibit 
the flagellar rotation rate [5] (thus effectively coupling the signal transduction system to 
the actuation system), though the precise mechanism through which this is achieved is 
unknown. An additional uncertainty lies in the demethylation connectivity between the 
CheB proteins and the two receptor clusters, though these questions have been the subjects 
of several studies, [5J HI [S]. Although how the CheB and CheY proteins interact with the 
sensing and actuation modules is not known for certain, it will be shown later in the paper 
that under some mild assumptions, FCD can be exhibited by the bacterium regardless of 
the exact structure of these connectivities. 

1.2 An MWC model of receptor dynamics 

We employ an MWC-type allosteric model for the receptor activities [T7]. Such models 
have been proposed for several bacterial chemotactic systems and have been found to be 
consistent with experimental data [T31 HH El [TU]. The main assumptions of the model 
are that receptors are either active or inactive, and that ligands have a higher affinity for 
inactive receptors than for active receptors. Respectively, wc denote by a(t) and a(t) the 
probabilities at time t of a transmembrane and cytoplasmic receptor being active. For each 
receptor, this probability can be approximated by the ratio of the Boltzmann factor of the 
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active state to the sum of the Boltzmann factors of all the states. Therefore if, at time i, 
the free-energy state of the membrane receptors is Ea when active and E\ when inactive, 
then the activity of the membrane receptors is approximated by 

a(t) = exp(-E A ) = 1 ,.s 

° U cxp(-E T ) +exp(-E A ) l + exp[-E A ) lj 

where Ea = Ej — Ea is the free energy difference between the active and inactive states. 

Similarly for the cytoplasmic receptors, the activity a(t) is dependent on their free-energy 
states when active and inactive, respectively Ea and Ei: 

a(t) = _ = L__ (5 ) 

exp(-Ei) + exp(-E A ) 1 + exp[-i?A] 

with Ea = Ei — Ea- The functions -Ea and Ea are assumed to have the same structure and 
take the form Ea = — [g m (m)+gL(L)] and Ea = — [gm{rh)+gL(L)] The functions g m ,gm are 
dependent on the methylation state of their respective receptors whilst the functions gh,gh 
quantify the effect of ligand binding on the receptor-free energy difference of the receptors. 
Following [IU] [7J [5] we make the assumption that each of g m ,g m is affinely dependent on 
the methylation state of its respective receptor cluster: 



n(m) = a(mo — m) and <7 m (m) = <S(mo — m) 



where a = a = 2 and mo = fho = 5. 

The binding of ligands to receptors leads to a loss of ligand translational entropy, propor- 
tional to the logarithm of the free ligand concentration [12l [10] . Due to the greater affinity 
of ligands to inactive receptors, this loss is greater in the case of ligands binding to active 
receptors. We denote the dissociation constants between ligands and active transmembrane 
(cytoplasmic) receptors by Ka (Ka), and between ligands and inactive transmembrane 
(cytoplasmic) receptors by Kj (Kj), with Ka 3> Ki and Ka 3> Ki due to the different 
affinities. From the E. coli chemotaxis literature, we adopt the values Kj = Kj = l8fiM, 
Ka = Ka = 3m M. As in [12], the change in receptor free energies due to ligand bind- 
ing to active transmembrane and cytoplasmic receptors is then, respectively, — ln(-^) and 

— ln(-~^). On the other hand the change in receptor free energy due to ligand binding to 
inactive transmembrane and cytoplasmic receptors is, respectively, — ln(^-) and — ln(-jp). 

The effect of this on the free energy differences Ea, Ea between active and inactive receptors 
can be characterized, as in [12] [TH], as 



9l(l) = to (i + i) - 1» (i + md 

fc( i )= to( 1+ !)-to(i + -|) 

for each cluster respectively. Due to the differences in affinities, we note that gh{L) and 
gL{L) are increasing functions of L and L respectively, which means that a(t) and a(t) are 
decreasing functions of L and L respectively. The greater affinity of ligands for inactive 
receptors therefore has the effect of shifting the receptors towards the inactive state. 
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Note that in the ligand concentration range Kj <C L <C K a and Kj < L < Ka, the 
receptor activities can be approximated by 



exp(a[m -m])-gj 
1 



and 



(6) 



exp(a[m - m])j^ 



2 Main results 

Following similar definitions in the literature, [SJ [T] we give the following definition of fold- 
change detection for the R. sphaeroides chemotaxis pathway. 

Definition 1. The R. sphaeroides chemotaxis system (QJ) exhibits fold change detection 
(FCD) in response to a sensed ligand input signal L(t) if its receptor activities a(t),a(t), 
initially at a steady state corresponding to L(0), are independent of linear scalings p > of 
the input L(t). 

Note that the chemotaxis protein phosphorylation network ([2]) takes as its sole inputs 
the signals a and a. For this reason, Definition [T] implies that if the system §S§ exhibits FCD 
in its activities, it also exhibits FCD in the concentration of its phosphorylatcd chemotaxis 
proteins (the elements of the vector x). The bacterium's flagellar behavior would also be 
expected to exhibit FCD as the flagcllum rotation rate is a function of the phosphorylatcd 
ChcY3, CI1CY4 and CheYg concentrations. Before giving the main result, we make the 
following assumption on the chemotaxis system dynamics. 

Assumption 2. The system has a unique steady state for any given L. 

Theorem 1. Under Assumptions [7] and and under approximation {6p the chemotaxis 
system with steady state initial conditions, will exhibit FCD in its activities a, a for 
ligand inputs in the range Ki <C L <C Ka cind Kj <C L <C Ka, in the sense of Definition]^ 

Proof. In the following, we assume that all ligand concentrations lie in the ranges Kj <C 
L <C Ka and Kj «I« Ka, and therefore approximation (|6]) holds. The proof is based 
on the existence of equivariances [I]. 

Suppose that in response to an external ligand input signal L = L\{t), the system (|3j). 
initially at a steady state corresponding to L = Li(0), exhibits a solution 
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and outputs oi(t) = G{m 1 {t), Lx{t)), L x (t) = C£i + D(u(xi, ai))L\, o x = G{rhi{t), L x (t)). 
Now if the ligand input is scaled to L = L2(t) = pLi(t), where p > 0, and if the initial state 
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is a solution of ([3]) since, under approximation ([5]), the outputs arc then 

a 2 = G{m 2 ,L 2 ) 

= G(mi + i logj5,pii) = G(mi,Li) = ai 

L 2 =C&+D{u{* 2 ,a 2 ))L 2 ' 

= P V C^ + D{u{yL U ar))p v L\ = p v L x 



a 2 =G{m 2 ,L 2 ) 

= G(fh 1 + ilogp v ,p v L 1 ) 



G(mi,Li) = ai 



which means that 



m 2 




rh 2 


d 






.6 . 





d_ rh 2 _d_ rhi(t) + j^fogp" 
dt x 2 dt xi(t) 

F(ai, w(xi)) 
#(ai,?i;(xi)) 
7J(xi,ai,ai) 
^[Mi+B(u(x.i, ai ))L^ 

F(a 2l w(-x 2 )) 
F(a 2 ,w(yi 2 )) 
H(yL 2 ,a 2 ,a 2 ) 

which verifies the claim that ([7]) is a solution of ^ when the ligand input is L{t) = pLi(i). 
Since the scaled inputs L = L\{t) and L = L 2 (t) = pLi^ yield the respective output 
pairs a\,a\ and a 2 ,a 2 and since a\ = a 2 and a\ = a 2 , it follows that system ([3]) under 
Assumption [1] exhibits fold change detection if the initial conditions of the system are rxii(O) 
when L = L\ and 1112(0) when L = L 2 (t) = pLi(t). In the language of [T], we have proved 
that the mapping mi i— > m 2 is an equivariance associated to scalar symmetries on inputs. 

Now if the system has a unique fixed point for any given L, and if mi(0) is the fixed point 
when L = Li(0), then ni2(0) is the fixed point when L = L 2 (0) = pLi(0) since, if rhi = 
when L = Li(0) then rii2 = when L = L 2 (0) = pL\(Q). Therefore if system ([3]) has a 
unique fixed point for any given L, starts from steady state conditions and yields solution 
mi(i), then scaling L by p > and initiating the system from steady state conditions will 
cause the system to yield the solution m 2 (t) and thereby exhibit FCD. □ 



3 Two R. sphaeroides chemotaxis models 

There are several integrated R. sphaeroides chemotaxis pathway models in the literature 
[UIHJIS]- In this section, we present two new models, differing from the previous ones in that 
their receptor dynamics arc of the form ((3]) and satisfy the MWC model given in Section 
11.21 Both of the new models we present were fitted to experimental data available in [J] and 
were able to reproduce the gene deletion data in [8l 0] . 

Each of the models presented satisfies the assumptions of Scction [TT2l and thereby exhibits 
FCD in the ligand range Ki <C L <C Ka and Kj <C L <C Ka- The dcmethylating feedback 
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structure for the models is restricted to that in 0] , which proposed an asymmetric feedback 
structure wherein CheBi demethylates both clusters and CheB2 demethylates the cytoplas- 
mic cluster, although a model with any feedback connectivity is capable of exhibiting FCD 
under the assumptions we make. 

The structural differences between the models lie in the signal L, which captures how 
external ligands are transduced to the cytoplasmic cluster. These models illustrate the point 
that, despite the differences in their internal connectivities, FCD behavior is conserved under 
the assumptions above. 

Following [Tni 0] , we make the assumptions that CheB proteins demethylate active re- 
ceptors, whilst CheR proteins methylate inactive receptors, and that CheR proteins operate 
at saturation. The CheR2 and CheR3 protein concentrations are therefore assumed to be 
constant and normalized to 1/iM each. Denoting by R 2 , R3 the concentrations of CheR2 and 
CheR3 and by Bi p ,B 2p the concentrations of phosphorylated chemotaxis proteins CheBi, 
ChcB 2 , mass action kinetics give the following general form for F, F in ([3]) 

m = F(a,a,w(x)) = k R (l - a)R 2 - k Bl B lp a - k B2 B 2p a 
m = F(a,a,w(ii)) = k R (l — a)R 3 — k B2 B 2p a 

where k R , k R > are methylation and k Bl , k B21 k B2 > demethylation rate constants. The 
probabilities of activity a, a given by ((H), ([5]). The models were obtained by fitting the 
constants kn,kn,k Bl ,k B2 ,k B2 in (JSJ). 

The experimentally measured output which was used to fit the model is the flagellar 
rotation frequency /. As shown in Figure [2] the CheY proteins control the rotation of the 
flagellum, and this is believed to happen through inhibitory binding [8|. The measured 
rotation frequencies to which we fit our models varied between Hz and a maximum of 
approximately 8 Hz. As discussed in [4], this maximum was very rarely exceeded, and is 
therefore assumed to be a physical limit on how fast the flagellum can rotate. As such, the 
rotation frequency is modeled as the Hill function 

/ 



0.125 + *(F s ,,y 4 „l'6,) 4 
where 

</>(Y 3p> Y 4p ,Y 6p ) = 0M2Y« 



' 0.1 + Y 3 +Ya 



(the negative sign denotes anti-clockwise rotation). In this way, / varies between 0-8 Hz, 
and decreases with increased concentrations of phosphorylated ChcY proteins. 



3.0.1 Model I 

Model structure: In this model the cytoplasmic receptors are assumed to sense internal- 
ized ligands, the concentrations of which are dependent on the external ligand concentration 
L. At the same time, as in [4], we assume there to be some interaction between the chemo- 
taxis proteins CheY3, CheY4 and the cytoplasmic cluster, and the function gh{L) takes as 
its input L = 10+ y 0L |_y- ■ ^ schematic of this model is shown in Figure |U 

A simulation of the model together with the tethered cell trace to which the model was 
fitted is shown in Figure [5l For comparison, Figure [5] additionally shows a simulation (with 
the same ligand input) of the model suggested in [4], which was fitted to the same tethered 
cell assay. The root mean squared error between the output of Model I and the tethered 
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cell assay is 0.88, which compares favorably to the corresponding error for the model in gj, 
which is 1.27. 

This model would be expected to exhibit FCD in the ligand ranges Ki <C L <C Ka and 
Ki ClC Ka- The latter range is equivalent to 

Kj (1 + 0.1[Y 3p + Y ip \) « L « K A (l + 0.1 [5^ + Y ip }) 

and since the total amounts of intracellular CheY3 and CheY4 (phosphorylated and un- 
phosphorylated) are 3.2/j.M and 13.2/jM respectively, then according to this model, simula- 
tions should show FCD in the range 2.64Kj <C L <C Ka- Figure [5] shows that this is indeed 
the case, with similar output traces obtained for the step changes in L from L = 1000/iM 
to 200 /iM and from L = 500/i M to lOO^M. 

Model parameters: k R = k R = 0.0045 k Bl = k B . 2 = 2.116 k B2 = 2.822. 



L 




Figure 4: Schematic of Model I. 



3.0.2 Model II 

Model structure: Here, the model's internally sensed ligands L are related to L via the 
differential equation L = —jL + jL. Whilst the ligand concentrations L, L modify receptor 
activities, the cytoplasmic receptors are otherwise unregulated by internal cell signals, and 
therefore L is not a function of u. A schematic is illustrated in Figure [71 and a simulation 
of the model together with the tethered cell trace to which the model was fitted is shown in 
Figure [8] For comparison, Figure [8] additionally shows a simulation (with the same ligand 
input) of the model suggested in [4], which was fitted to the same tethered cell assay. The 
root mean squared error between the output of Model II and the tethered cell assay is 0.95, 
which, as with Model I, also compares favorably to the corresponding error for the model 
in gj, which is 1.27. 

Note that if L were to undergo a step change from L = L a fiM to L = LbfiM and if the 
system is initially at steady state (where L(0) = L a ), then L would remain confined to the 
set [L a , Lb). Therefore, for such a step change, Kj <L< Ka implies that Kj < L C Ka- 
Figure [S] shows that simulations of this model do show FCD in this input range, with similar 
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Figure 5: Simulation of Model I (red) in response to a step rise (at 245 seconds) and fall (at 
370 seconds) in the ligand level L from L = to L = 100 and back to L = 0, with a tethered 
cell assay (black). The dashed blue trace is a simulation of the previously published model 
in [1] subject to the same ligand input. 



output traces obtained for the step changes in L from L = 1000/iM to 200 fjM and from 
L = 500/i M to 100^M. 

Model parameters: k R = k R = 0.0057 k Bl = k B . 2 = 2.376 k B2 = 2.970. 
3.1 Future experiments for model invalidation 

The models presented above are two systems based on the assumptions of Section 11.21 that 
reproduce the experimental data of [3], but which additionally show FCD. By comparison, 
the model suggested in [4], based on different receptor dynamics, does not exhibit FCD in 
response to the inputs used in the simulation in Figures [6j [9] as shown in Figure [lOj 

One important feature to note is that the FCD property is preserved regardless of the 
exact dynamics in ([2]) and regardless of the interactions between the receptors and the 
chemotaxis proteins, as long as the conditions of Theorem [T] arc satisfied. Therefore, if the 
receptor dynamics model given in Section 11.21 is accurate, then in the ligand concentration 
ranges Ki <C L <C Ka and Kj « L « Ka, FCD is a robust dynamic property of the 
chemotaxis system that should be observed in both wild type and in mutant strains of R. 
sphaeroides that have chemotaxis protein deletions and over-expressions. 

The above points suggest experiments in which the FCD dynamic phenotype can be 
used to discriminate between Models I and II on the one hand, and the model suggested in 
[4] on the other: 

• If Models I and II are to invalidate that of [3] then the wild type bacterium, initially at 
steady state, should show near identical flagellar output responses to the step ligand 
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Figure 6: Model I output in response to step changes in L from L — 1000/iM to 200 fjM 
and from L = 500^ M to 100/xM 




inputs L = 1000^M to 200 fiM and L = 500/x M to 100. 

• Overcxprcssing the chemotaxis protein CheY4 five fold was shown in [8] to not destroy 
the chcmotactic response of the bacterium. Such a mutant strain should therefore, 
according to Models I and II, also exhibit FCD in response to a range of step changes 
in the external ligand concentration L. We can calculate this range for each of the 
two models as in Sections 13.0.11 and 13.0.21 In Model I, the five-fold increase in CheY4 
means that FCD should be observed within the range 7.92Kj <C L <C Ka, whereas 
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Figure 8: Simulation of Model II (red) in response to a step rise (at 245 seconds) and fall (at 
370 seconds) in the ligand level L from L = to L = 100 and back to L = 0, with a tethered 
cell assay (black). The dashed blue trace is a simulation of the previously published model 
in [1] subject to the same ligand input. 



for Model II this range is Kj < L < Ka- 

• If we define the adaptation time when the model is subject to a step decrease in 
ligand to be the time that it takes from the application of the step for the deviation 
of the flagellar rotation frequency from its steady-state value to fall to 25% of its 
maximum, then under this definition, the adaptation times for Models I and II are 
65 seconds and 62 seconds respectively, whereas for the model in [JJ, the adaptation 
times are 266 seconds for the step ligand concentration decrease of 1000 fj, M to 200 
\x M, and 162 seconds for the step decrease of 500 /x M to 100 fjM. If Models I and 
II are to invalidate those in [4], the experiments should yield approximately equal 
adaptation times in response to these two step changes in ligand concentration, and 
these adaptation times should be around 60 seconds. 

Further model discrimination between Models I and II can be performed using the tools 
presented in [5J H] . 

4 Discussion 

The models presented herein differ from earlier R. sphaeroides chemotaxis models in two 
main respects: first, the receptor dynamics are based on the MWC allosteric model. This 
model has been shown to be a fairly accurate representation of the receptor dynamics in E. 
coli. The homologies between the bacteria and the similarities between their overall chemo- 
taxis mechanisms give us reason to believe that the MWC model may, under experimental 
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Figure 9: Model II output in response to step changes in L from L = 1000/iM to 200 fjM 
and from L = 500^ M to 100/dVI 
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Figure 10: Simulations of the model in [4], subject to step changes in L from L = 1000/^M 
to 200 fiM and from L = 500^ M to 100/dVl 



testing, eventually prove to be a realistic way of representing the R. sphaeroides receptor 
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dynamics. 

The second point of departure of these models from earlier ones is that the assumptions 
on the possible relationships between the external and internal ligand concentrations are 
relaxed to admit dynamic relations. The motivation behind this model is to capture any 
phase delays between sensed changes in the external ligand concentration and the effect of 
such changes on the internal cell environment. 

The external-internal ligand relation of Model I closely follows that of [4]. In effect, the 
activity of the cytoplasmic cluster depends on the external ligand concentration, L, and, 
indirectly, on the activity a of the membrane cluster via the phosphorylated chemotaxis 
proteins CheY3-P and CheY4-P, as schematically illustrated in Figure 2] On the other 
hand, the cytoplasmic receptor activity in Model II does not depend on any chemotaxis 
proteins, and its sensed ligand signals are merely phase-delayed versions of the external 
ligand concentration. 

As experimentally shown [5], chemotaxis requires CheYg and one of either CheY3 or 
CI1CY4, as deletion of either CheYg or both of CheY3 and CheY4 destroys the chemotac- 
tic ability of the bacterium. In the models we present, this was captured by the interaction 
of the three CheY proteins at the flagellum in what is effectively an AND logic gate that 
will only be activated if both CheYg and at least one of CheY3 p or CheY4 p are present. The 
signal transduction dynamics (Q31IE]) snow that CheY 3p and CheY4 p are solely phospho- 
rylated by the membrane cluster, whereas CheYg receives most of its phosphates from the 
cytoplasmic cluster. In essence, this structure means that there are essentially two paths 
from the external ligands to the flagellum that terminate at the AND gate: one path via 
the membrane cluster in which CheY3 and CheY4 proteins convey the signal, and one 
path via the cytoplasmic cluster, in which CheY6 p conveys the signal. This resembles a 
recurring biochemical motif [11], and the selective advantage it bestows could be improved 
energy taxis [3] with respect to simpler chemotaxis circuits such as that of E. coli. The main 
feature of this improved pathway is that the flagellar motion will only vary if both signalling 
paths from L to the flagellum are activated. Since the cytoplasmic cluster may integrate 
un-modeled metabolic information from within the cell, it would be important that any 
variation in flagellar activity only results from a change in the metabolic state of the cell 
that arises from a change in the local chemoeffector environment. If this is indeed the case 
then the signalling path from the cytoplasmic cluster is only activated if the metabolic state 
of the cell changes, whilst the signalling path from the membrane cluster is only act.imt.otl 
if the immediate chemical environment changes. Only if both are activated together would 
the cell 'know' that the change in its metabolic state is due to a change in chemoeffector 
concentration, and only then would it change its flagellar activity. 

4.1 Selective advantage of FCD 

Whether FCD bestows upon the bacterium a selective advantage or simply arises as a 
by-product of the chemotaxis system's structure is a question of interest. The fact that 
FCD is present in simpler chemotaxis circuits than that of R. sphaeroides (e.g. in E. coli) 
suggests that the advantages gained by having such a property would be independent of 
the complexity of the bacterium's chemotaxis pathway. It may be that the metabolic payoff 
to the bacterium of moving to more chemically favorable regions depends on the relative 
chemical improvement in its environment rather than the absolute change. A reason for 
this could be that biasing its movement towards longer swims could be metabolically costly 
for the bacterium, and moving in this way is only worthwhile if the metabolic gain is 
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significant. A potential disadvantage of FCD to the bacterium could be a high sensitivity 
to small fluctuations in sensed ligand when the background ligand concentration is low, due 
to the fact that the gain in the flagellar rotation frequency would then be high. However, 
this disadvantage is offset by the fact that FCD behavior only occurs at background ligand 
concentrations significantly above a threshold, given by Kj in the models above. 

4.2 FCD as dynamic phenotype for model invalidation 

The models we have presented provide an example of how dynamic phenotypes can be used 
to discriminate between competing biochemical models. Given two models of the same sys- 
tem, a mathematical analysis can be used to identify regions in the parameter and input 
spaces in which a certain qualitative dynamic behavior, such as FCD, could be expected. 
Ideally, this behavior would be expected to be robust to any genetic mutations or envi- 
ronmental conditions, and the conditions under which this behavior would occur would be 
implcmcntable experimentally. Model discrimination can then be performed on the basis of 
whether or not the system robustly reproduces the dynamic phenotype experimentally. This 
differs from traditional forms of model discrimination in that it can be used to discriminate 
between different biological mechanisms, and can be used to identify whether an observed 
phenomenon is due to the fine tuning of biological parameters or due to a more fundamental 
structural property of the system. 
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